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I.  Introduction 

The  unalyueal  dowhward  continuation  (ADC)  is  widely  used  in  the  solving  of  the  geodetic 
boundary  value  problem.  One  example  is  the  use  of  the  ADC  for  the  solving  of  Molodensky  s 
problem  (Moritz,  1980,  Section  45).  Recently,  the  airborne  and  satellite  gradiometry  have  been 
getting  substantial  improvement  and  they  will  supply  a  huge  gravity  gradient  data  for  the 
determination  of  the  Earth's  gravity  field.  To  process  such  data  the  ADC  method  has  to  be  used 
and  sonic  lefercnce  surface,  such  as  the  sea  level,  the  ellipsoid,  may  be  chosen.  The  topographic- 
effect  is  becoming  an  important  problem  in  the  use  of  the  gradient  data.  This  study  is  addressed  on 
the  problem  of  the  topographic  effect  on  the  using  of  the  ADC  method. 

In  Chapter  2  the  classic  topographic  reductions  -  the  Bouguer  reduction  and  the  Helmert 
second  condensation  will  be  introduced  into  the  airborne  gradiometry.  The  reductions  are  studied 
in  the  space  domain  and  the  frequency  domain.  The  reduction  in  space  domain  is  tedious  but  they 
g’vc  a  clear  view  of  the  physical  meaning.  In  the  frequency  domain,  this  problem  becomes  much 
easier  to  solve. 

In  Chapter  3  die  topographic  effect  for  the  use  of  the  analytical  downward  continuation 
method  is  studied.  The  remove-restore  of  the  mass  above  the  ellipsoid  is  used.  The  goal  of  the 
study  is  to  find  the  topographic  corrections,  in  another  word,  to  find  the  difference  between  the 
disturbing  potential  on  the  sea  level  and  disturbing  potential  which  is  downward  continued  to  the 
sea  level  from  the  outside  of  the  earth. 

In  Chapter  4  we  study  the  topographic  effect  on  the  solving  of  the  Molodensky's  problem.  It 
was  proved  that  the  ADC  solution  is  equivalent  to  the  Molodensky's  solution  <ibid,  p.  388).  In 
this  chapter  we  show  that  the  topographic  effect  still  exists  in  the  solution  of  Molodensky's 
problem  by  using  the  ADC  method,  even  if  this  effect  is  very  small.  The  topographic  effect  on  the 
determination  of  the  coefficients  of  the  spherical  harmonics  of  the  geopotential  by  using  the  ADC 
method  is  also  considered. 


2.  Topographic  Correction  on  the  (Iravity  Gradient  Data  Measured  by  an 

Airborne  Gradiometer 

2.1  Introduction 

In  this  chapter  we  consider  the  topographic  correction  on  the  gravity  gradient  data  measured 
by  an  airborne  gradiometer.  W'e  assume  that  the  airborne  flights  are  at  a  constant  altitude  above  the 
sea  level.  In  comparison  with  the  satellite  gradiometry  the  airborne  flights  are  much  lower  and  can 
sense  much  higher  frequencies  of  the  gravity  field.  It  is  well  known  that  the  gravity  gradient 
depends  strongly  on  the  roughness  of  the  topography.  Because  the  airborne  gradiometer  flies 
close  to  the  topography,  we  will  try  to  eliminate  the  topography  effect  directly  from  the  measured 
gravity  gradient  data. 

The  topographic  correction  for  airborne  gravity  and  gradient  data  has  been  discussed 
extensively  in  the  literature;  e.g.,  Chinnery  (1961),  Parker  (1972),  Dorman  and  Lewis  (1974), 
Hammer  (lQ7fi)  Recently  an  algorithm  was  developed  for  eliminating  the  "topographic  noise" 
(Tziavos  et  al.,  1988).  The  fast  Fourier  transformation  was  used  and  this  algorithm  is  very 
efficient  to  compute  the  topographic  effect  on  the  gravity  and  gradient  data. 

Now  we  want  to  make  the  statement  "eliminating  the  effect  of  the  topography"  more  clear. 
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Tziavos  et  al.,  (1988)  developed  an  algorithm  to  compute  the  gradient  of  the  attraction  of  the 
topography.  The  topographic  effect  was  considered  as  "topographic  noise"  (ibid)  and  subtracted 
from  the  measured  gradient  data. 

Actually  the  topographic  effect  is  a  part  of  the  gravitational  potential  of  the  earth  and  cannot  be 
simply  eliminated  by  subtracting  its  contribution  from  the  measured  gravity  and  gradient  data. 
Before  we  consider  the  effect  of  the  eanh's  topography,  we  must  ask  ourselves,  what  is  the  goal 
and  what  will  be  determined  by  using  the  gravity  gradient  data  or  the  reduced  data.  This  is  very 
important  for  the  start  of  the  topography  reduction  procedure. 

If  the  contribution  of  the  topography  is  subtracted  from  the  gravity  and  gradient  data,  it  means 
the  topography  mass  is  removed,  then  the  disturbing  potential  and  its  derivatives  are  harmonic 
above  the  sea  level  and  can  be  analytically  downward  continued  to  the  sea  level  without  any 
theoretical  difficulty.  Obviously,  the  reduced  gravity  anomaly  is  the  Bouguer  anomaly  and  the 
topography  reduction  is  the  Bouguer  reduction. 

In  this  chapter,  we  assume  that  the  gcoid  and  the  deflections  of  the  vertical  will  be  determined 
in  a  mountain  area  by  using  the  aerial  gravity  gradient  data.  The  topographic  effect  will  be 
eliminated  from  such  data.  For  this  purpose  an  instinctive  action  is  to  expand  the  Bouguer 
reduction  and  the  Helmert's  second  condensation  into  the  airborne  gradiometry. 

2.2  Topographic  Correction  Studied  in  Space  Domain 

2.2.1  Bouguer  Reduction 

The  Bouguer  reduction  in  airborne  gradiometry  is  to  eliminate  the  "topographic  noise"  by 
subtracting  the  contribution  of  the  topography  to  the  gravity  and  gradient  data  (e.g.,  Tziavos  et  al.. 
1988).  The  airborne  ^diometry  is  taken  in  a  local  area.  For  such  area  a  plane  approximation  for 
the  computation  is  suitable  (Jekeli  et  al..  1985). 

The  geometry  of  the  computation  of  the  topographic  effect  in  the  airborne  gravity  gradiometry 
is  illustrated  in  Fig.  1. 
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Fig.  1.  Geometry  of  the  Computation  of  the  Topographic  Correction 
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We  define  the  lixal  coordinate  system  whose  origin  is  located  at  sea  level  and  the  x  axis  points 
to  the  north  and  y  axis  points  to  the  east. 

The  potential  of  the  mass  above  the  sea  level  at  the  point  Q  is  given  by: 


with 

L"=  s‘'  +  (zo- z)^ 


(2-n 


(2-2) 


where  G  is  the  gravitational  constant  and  the  density  of  the  topography  p  (x,  y.  z)  is  a  function  of 
the  location  of  the  current  point,  x  is  a  2D  plane  which  coincides  with  the  local  sea  level  and  s-  = 
(x-xq)“  -1-  (y-yQ)2,  At  first  we  assume  that  the  density  p  is  a  constant.  Notice  that 


s‘  -1-  Zn>  2zn  z  -  z^ 


(2-3) 


then  eq.  (2- 1 )  can  be  expanded  by 

Vq=Gp  ff  ^dxdy^icp  // 


.3  rr  r‘(2zoz-z ) 

-  (2-4) 

where  L'  ^  +•  z^,  is  the  distance  between  the  point  Q  and  the  current  point  on  the  sea  level  x. 

The  equation  (2-4)  gives  the  potential  of  the  toposcaph)  at  th?  point  Q.  Th**  'erms 
contain  the  inverse  of  the  power  of  L'  higher  than  5  are  omitted.  In  the  next  discussion  we 
separate  the  topography  into  a  Bouguer  plate  with  thickness  hp  and  a  terrain  undulation  around  the 
point  P. 

TTe  Bouguer  plate  generates  a  potential  at  the  point  Q: 

Vq=Gp//^  //rdxdy 


/'*  2Za  Z'i. 

- dzdxdy  + 

5 


'3 


=  ImJtGphpS  -  2rf5p(;„hp-  ihj) 

where  denotes  the  potential  of  the  Bouguer  plate.  The  first  term  in  (2-5)  becomes  infinite  when 
s— >oo.  Tnis  infinity  is  due  to  the  plane  approximation.  In  spherical  approximation  this  term  is  a 
constant  (cf.  Moritz,  1968).  From  (2-5)  we  can  see  that  the  potential  of  the  Bouguer  plate  is  a 
function  of  the  flight  altitude,  zp.  Obviously,  the  potential  of  the  Bouguer  plate  is  not  a  regular 
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potential  because  it  does  not  s  anish  at  infinity.  The  terrain  undulation  has  the  contribution  to  the 
geopsrtential  at  the  point  Q: 


dxdv 


=  Gp  fj  dxdy  +  I  Gp  fj  J  '  dxdy 


'T* 

where  Vq  denotes  the  potential  of  the  terrain  undulation  around  ih**  point  Q. 

Now  we  consider  the  first  and  second  derivatives  of  the  potential  of  '^e  topography. 
The  senical  derivative  of  the  potential  of  the  topography  is  given  by 

v.=  !f's=Gp///^  ^  ' 

azo  ^  •'o 


0  Biq 


L 


dzdxdy 


=  V  +  V 
’  2  ’  2 


with 


V,  =  Gp  f  f  f  '■ — -pdzdxd 
•'  •'t  *^0  dZo  ^ 

v;=Gp///‘Ai 

J  •'hp  az-  ^ 


dzdxdv 


p  ozo 


(2-6) 


(2-7) 


(2-8i 


(2-9) 


Obviously.  V^  is  the  contribution  of  the  Bouguer  plate  and  is  the  contribution  of  the  terrain 
undulation  Notice  that 


//7= 
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“  dxdv 


V  s“+|z„-a) 

=  27t  lim  s  -  27t  (z.  -  a') 
s— >  c»  V  u  / 

where  a  >  0  is  a  parameter,  then  we  obtain 


"  V s’Tii-a) 


(2-10) 
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=^Gp 


1 


dxdy 


--27iGphp  ^2-11) 

This  is  nothing  else  than  the  attraction  of  the  Bouguer  plate,  where  the  computation  point  is  not 
only  on  the  Bouguer  plate  but  also  at  the  flight  altitude.  Equation  (2-11)  means  that  the  attraction 
of  a  Bouguer  plate  is  a  constant  on  and  above  the  Bouguer  plate  and  it  does  not  change  with  the 
flight  altitude  zp. 

The  effect  of  the  terrain  undulation  is  given  by 


\ 


h^-hj 

- dxdv 


'3 


(2-121 


The  attraction  of  the  mass  above  the  sea  level  at  the  point  0  is  then 

/  ->\ 


V. 


2:tGphp 


■  V 

I, 


(2-13i 


In  a  mountainous  area  the  first  term  in  (213)  is  the  principal  value  and  much  bigger  than  the  last 
term  in  eq.  (2- 1 3). 

In  the  same  way  we  can  get  the  honzontal  derivatives  of  the  potential  of  the  topography  The 
horizontal  denvatives  of  the  potential  of  the  topography  are 


—  ^ d/dxdv 


- T-  dzdxdv 


We  define  the  horizontal  derivatives  of  the  potential  of  the  Bouguer  plate  as; 


(2-14) 


(2-15) 


(2- 16a) 
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SVp  ^  J  n  h  (^2-’ 

It  IS  etisv  to  prove  that  the  Bouguer  plate  has  no  effect  on  the  horizontal  derivatives; 

Introduce  the  polar  crxirdinate  system; 

X  =  X  -  Xp  -  s  cos  6 

'i'  =  \  -  Np  =  s  sin  6  1 

so  ihat  has'e 
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’’dxdv 
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=  — (-27iGphp)  =  0 

dzp  (2-22) 

Based  on  the  above  discussion  we  can  see  that  the  Bouguer  plate  has  no  effect  on  the  gradient 
data.  Therefore,  it  is  not  possible  to  recover  the  effect  of  the  Bouguer  plate  on  gravity  and  the 
potential  by  processing  the  aerial  gravity  gradient  data  only. 

By  using  eqs.  (2-18)  and  (2-22)  we  have  the  second  derivatives  of  the  potential  of  the 
topography 
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Q  _ 


axpazo 

aV 


=  -3Gp 


x-x 


Q  _ 


3yp3zo 


=  -3Gp 


fu: 


ff 


fJM 


pK^o-z) 


dzdxdv 


dzdxdy 


T  "P 


(2-26) 


(2-27) 


v„-4--op//r 


dz 


3z 


dxdy  +  3  Gp 


ffP- 

J  P  •'hp  L 


dzdxdy 


(2-28) 


Between  the  six  components  V^y . there  are  five  independent  components.  A  strict 

condition  is 


Vxx  +  Vyy  +  =  0 


(2-29) 
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The  derivation  of  the  formulas  are  routine  and  tedious,  but  it  helps  us  to  understand  better  the 
physical  meaning  of  the  reduction  and  helps  us  to  figure  out  the  role  of  the  Bouguer  plate  and  the 
terrain  undulation  in  the  airborne  gravity  gradiometry. 

Next  we  consider  a  more  general  case  where  the  density  of  the  topography  is  not  a  constant 
We  separate  the  density  function  into  a  constant  part  and  a  variable  part: 


p(x,y,z)  =  Po+8p  (x,y.zl 


(2-30) 


where  po  is  a  constant  and  5p  is  the  variable  part  of  the  density  of  the  topography  which  is 
sometimes  called  the  density  anomaly  (Forsberg,  1984).  We  assume  that  this  density  function  is 
known.  For  the  constant  density  pg  all  formulas  derived  above  remain  valid.  Therefore  we  need 
only  to  consider  the  potential  caused  by  the  variable  part  of  the  density. 

The  potential  caused  by  -c  density  anomaly  5p  is  given  by 

h8p 


L 


dzdxdv 


(2-31) 


AV, 


1 


1  2zoZ-z 


'2 


dzdxdy 


Define  the  new  functions 
.  h 


5p^,=  I  5p(x,y,z|dz 
•'o 


5p 


I  (2zgz-  z^)sp(x.  y,  z| 

•^0 


dz 


Then  (2-32)  can  be  written  as 
1 


AVq-  G  5p.; 


L' 


:-G6p^.. 


(2-32) 


(2-33) 


(2-34) 


(2-35) 


where  "  denotes  the  convolution  of  two  functions.  Therefore  the  numerical  evaluation  of  (2-35) 
can  be  done  by  using  the  efficient  fast  Fourier  transformation.  The  first  derivatives  of  the  potential 
AVq  are  given  by: 


dAVo 
AV .  = - 2 

ax„ 


-// 


x-x. 


L 


^  5pM(’‘>y)axdy-t-2  G 


If. 


x-x. 


5pn(x,  yjdxdy 


(2-36) 
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^\’  = 


9AV 

3yp 

aAVo 

AV  = - ^=-Gz, 

az. 


^  ^  J  J  “7^  Spfjjx,  y)  dxdy 

L  I-* 


of  f  dxdy  -  ^  G  zo  f  j  dxdy 


(2-37) 


(2-38) 


If  the  density  is  not  a  constant,  then  the  second  derivatives  of  the  Bouguer  plate  is  no  longer 
equal  to  zero.  Therefore  we  have; 


(x-XpKy-yp| 


SpM  dxdy 


(2-39) 


AV 


xz 


2 

a  AV 


ax  ay 


S=-3Gzo 


X-Xn 

-^5P|^  dxdy- 
L  ^ 


x-x. 


•  5p^  dxdy 


(2-40) 


(2-41) 


(2-42) 
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a  Av 

-  avpazp 


-  3G 


J^off  “Sp^dxdy 

•’x  L 


(2-43) 


(2-44) 


The  inverse  of  the  higher  powers  of  the  distance  L'  appearing  in  the  formulas  means  that  the 
derivatives,  such  as  AV^*,  AVyy,  ....  depend  strongly  on  the  value  of  the  innermost  zone.  The 
effect  of  the  remote  zone  is  very  smalt  and  vanishes  very  fast.  In  the  numerical  computation  the 
integration  region  can  be  taken  small,  but  the  data  should  be  dense  enough  around  the  computation 
point,  in  order  to  ensure  the  computation  accuracy. 

The  higher  power  of  the  inverse  of  the  distance  L'  are  omitted  in  the  above  formulas.  These 
higher  power  terms  contain  very  high  frequencies  of  the  potential.  In  the  numerical  computation 
the  omitted  terms  may  not  cause  serious  error  in  the  results. 

The  formulas  derived  above  are  in  convolution  form  and  the  efficient  numerical  computation 
method  -  fast  Fourier  transformation  can  be  used. 


2,2.2  Heltnert's  Second  Condensation 

It  is  a  well-known  fact  that  the  indirect  effect  of  the  Bouguer  reduction  on  the  geoid  is  very- 
large  (Heiskanen  and  Montz,  1967,  p.  142)  and  it  is  not  suitable  for  the  determination  of  the  geoid. 

In  this  paragraph  we  introduce  the  Helmert  s  .second  condensation  into  airborne  gradiometry 
reduction.  The  Helmen's  second  condensation  is  to  press  the  topographic  mass  onto  the  geoid,  or, 
in  a  good  approximation,  onto  the  sea  level. 

We  assume  the  density  is  constant,  so  that  the  density  of  the  condensed  layer  is  ph. 

The  potential  of  the  condensed  layer  at  the  point  Q  is  given  by  (cf.  Figure  1): 

Wq=G  ff~  dxdy 


First  we  consider  the  conden.sation  of  a  Bouguer  plate  with  thickness  hp  into  a  layer. 


(2-45) 
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The  density  of  the  layer  condensed  from  a  Bouguer  plate  is  p  hp  and  its  potential  at  the  point  Q 
is 


Wq=Gp  dxdy 

=  27cGph_lim  s-2jTz.h„ 

Ps-»oo  0  P  (2-46) 

Combining  (2-5)  and  (2-46)  we  obtain  the  variation  (change)  of  the  potential  of  the  Bouguer 
plate  due  to  Helmert's  condensation: 

B  R  R 

8Vq=  -  Vq+  Wq 

=  -JtGphp  ^2-47) 

R 

The  variation  of  the  potential  oVq  is  independent  of  the  flight  altitude  and  is  a  constant  above 
and  on  the  Bouguer  plate.  Eq.  (2-47)  means  that  if  a  Bouguer  plate  is  pressed  into  a  layer,  the 
potential  will  change  at  point  Q  a  constant  -  JiGphp.  This  term  is  small.  If  we  take  p  =  2.67 
g/cm3,  hp  =  1000  meters,  the  term  -  ttGphp  causes  a  geoid  undulation  change  of  5  cm. 

Roughly  speaking,  the  condensed  Bouguer  layer  generates  almost  the  same  potential  as  the 
Bouguer  plate  does.  The  difference  JiGphn  is  a  small  quantity  and  it  can  be  neglected  for  the  first 
approximation.  But  if  the  geoid  is  needecf  to  be  precisely  determined,  this  term  is  also  imponant 
and  cannot  be  omitted. 

Equation  (2-45)  can  be  written  as 

Wq=Gp  II  dxdv  -fGp  II  — E  dxdy 


(2-48) 


Combining  eqs.  (2-5),  (2-6)  and  (2-46),  (2-48)  we  get  the  change  (at  point  Q)  of  the  potential  of 
the  topography  after  Helmen's  condensation: 


8Vq=-Vq+W(; 


=  -rtGp  h_- 


lop/j; 


'3 


dxdy 


(2-49) 


The  last  term  in  (2-49)  is  one  order  smaller  than  the  first  term  and  can  be  neglected  except  in  some 
extreme  cases. 


Now  let  us  consider  the  variation  of  the  derivatives  of  the  potential  due  to  the  Helmen  s 
second  condensation. 

The  venical  derivative  of  the  potential  of  the  condensed  topography  at  point  Q  is 


=  ^^=-GpZo  f  f  — dxdy 
^’0  L  ^ 


We  separate  (2-50)  into 


W,  =  wf-t- wj 


(2-50) 


(2-51) 


'^--GpzoJJ  -^^dxdy 


wl=- 


=  -Gpzo// 


(2-5^) 


(2-53) 


Obviously  (2-52)  is  the  attraction  of  the  condensed  Bouguer  plate  and  (2-53)  is  the  attraction  of  the 
condensed  terrain  variations  with  respect  to  the  Bouguer  plate. 


.Nonce  that 


I  f  —  dxdy=  r^^dO  f 
•J  Jr  I  3  J  Jil 


(s^-t-zj) 


(2-54) 


we  have 


W,  =  -  2rtGph, 


(2-55) 


Comparing  (2-55)  with  (2- 1 1 ),  we  find  that  the  conden.sed  Bouguer  plate  layer  generates  the  same 
attraction  as  the  Bouguer  plate  does. 

Therefore  the  attraction  of  the  condensed  topography  is  given  by  combining  (2-51),  (2-53) 
and  (2-55); 


=  -  2TtGph  p  -  Gpzo  f  f  — dxdy  , 
J  Jj  I  '3 


(2-56) 
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and  the  difference  in  attraction  between  the  actual  potential  and  tne  potential  of  the  condensed 
topography  is 


1  r  r  h  -h; 

5V^=  -  V,+ W,  =  jG  pj  J  — dxdy 

^  L  ''  (2-57) 

It  is  interesting  to  compare  (2-57)  with  the  result  of  Vanicek  and  Kleusberg  (1987,  eq.  (14)).  They 
are  almost  the  same.  The  difference  is  just  due  to  the  location  of  the  computation  point.  If  the 
computation  point  Q  on  the  flight  altitude  is  down  to  the  point  P  on  the  topography,  then  eq.  (2-57) 
is  reduced  to  the  Vanicek-Kleusberg  result.  If  the  computation  point  Q  is  down  to  the  condensed 
layer,  is  reduced  to  the  classic  terrain  correction  (cf.  Wang  and  Rapp,  1989). 

Because  we  have: 


ff  ^dxdy=  r^"sinede  f” 

L  ^  -'o  -'o 


s"  ds 


11  2  2^ 

(s  -t-Zo) 


=  0. 


s^ds 
3/2 


=  0, 


(2-.58) 


(2-59) 


then  we  have  the  horizontal  derivatives  of  the  potential  of  the  condensed  topography  can  be  written 
as: 


9Wn 
W,  =  — 3 

9x„ 


L 


x-x 


pjdxdy 


9W( 

3yp 


=  Gp  f  f  (y-yp)dxdy 

^  3 


(2-60) 


(2-61) 


Combining  eqs.  (2-20)  and  (2-21)  with  (2-60)  and  (2-61)  we  get  the  change  (actual  topograph) 
minus  condensed)  of  the  horizontal  derivatives: 


8V,  =  -V,-.W,  =  -|gp  f f 


^o(h^-hp)4(h’-hJ) 


'5 


(x-Xp)dxdy 


(2-62) 


6Vy=-  Vy-HWy  = 


'5 


-(y-Ypjdxdy 


(2-63) 


lA 


The  vanations  are  one  order  smaller  than  the  horizontal  derivatives  of  the  potential  of  the 
topography.  It  is  easy  to  show  that  the  condensed  Bouguer  plate  has  no  effect  on  the  gradient  data. 

From  eqs.  (2-58)  and  (2-59)  we  know  that  the  second  derivatives  of  the  potential  of  a 
condensed  Bouguer  plate,  such  as  vV^r.  Wpy,  ....  defined  by 


(2-64a) 


(2-64b) 


are  equal  to  zero. 

D 

For  the  component  we  have 


2 


— 2.  dxdy 
L  ^ 


=  -  2kGp  h„—  (l)  =  0 


(2-65) 


Here  we  need  to  point  out  that  the  elevation  hp  is  considered  as  a  constant  in  all  integrals. 
Obviouslv,  hp  is  a  function  of  the  point  P,  therefore  the  order  of  the  elevation  hp  with  the  partial 
differential,  such  as  it  is  in  (2-64a)  and  (2-64b),  cannot  be  changed.  By  way  of  exception  the 
venical  differential  0/9zo  and  are  changable  with  the  hp  because  the  elevation  hp  (xp,  yp)  is 
not  a  function  of  the  flight  altitude  zq.  The  second  derivatives  of  the  potential  of  the  condensed 
topography  are  given  by 


W„  =  -^WQ=-Gpr  f  ^dxdy-^3Gpr  f  dxdy 

axp  l’  l' 


(2-66) 


2 

^Wq=3Gp  rr  ^|x-xpj(y-yp)dxdy 


(2-67) 
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W..,= 


yy  2 

3yo 


WQ=-Gpyj^  dxdy  +  3Gp J J  — jx-Xp)^  dxdy 


,2 

l' 


(2-68) 


w„= 


dXpdzQ 


Wq=-3Gpzo  J -5^(x-Xp)  dxdy 


(2-69) 


W  = 

’’  yz 


5xp9zo 


Wq=  -  3  Gp  zo  f  f  ^(y-yp)  dxdy 


(2-70) 


W,  = 


BZr 


Wq=-GpJ'^  dxdy  +  3GpJ'J^  -!^zodxdy 


h_h„  2 

L 


(2-71) 


By  combining  the  above  formulas  (2-66)  -  (2-71)  with  the  second  derivatives  of  the  potential  of  the 
topography,  eqs.  (2-23)  -  (2-28),  we  obtain  the  change  of  the  second  derivatives  of  the  potential 
due  to  the  topography  condensation: 


=  -V,,+  W,,  =  |Gprr  ii-dxdy-^Gprr  toLndxdy 

J  J  5  -  ‘Z  J ^  7 


yOpJ  ~(x-Xp)|y-yp)dxdy 


S^y=-^y 


^^yy  =  j^pjf^  -^dxdy-yGpJJ^  H  dxdy 


=  -jGpf  I  - ^(x-Xpjdxdv -t- 

“  -'t  L  5 


^  z5(h^-hp)-Zo(h^-hj)  +  ^(h'‘-hj) 


dxdy 


(2-72) 


(2-73) 


(2-74) 


6v  =  -  V  +w  =  - 

yz  *  yz  ^  yz 


h^-hj 

-'x-Xp)dxdy  -t- 


(2-75) 
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x-x. 


'7 


(h^-h5)-2o(h^-hJl 


dxdy 


w. 


(2-76) 


(2-77) 


(2-78) 


The  above  formulas  are  complicated.  But  they  can  be  written  in  the  convolution  form  and  the  fast 
computation  method  -  fast  Fourier  transformation  can  be  applied.  Now  wc  consider  the  case  in 
which  the  density  of  the  topography  is  not  a  constant. 

Following  eq.  (2-30)  we  separate  the  density  of  the  topography  into  a  constant  part  and  a 
variaiable  part.  We  consider  only  the  effect  of  the  variable  density.  Using  Helmert’s  second 
condensation  we  condense  the  topography  whose  density  varies  into  the  sea  level  with  the  surface 
density  5pM  defined  by  eq.  (2-33^  Then  the  potential  of  this  layer  at  point  Q  is 


^Wq  =  G  j  I  ^  d.xdy 

J  h  (2-79) 

Comparing  (2-79)  with  (2-32)  we  find  the  change  of  the  potential  due  to  the  condensation  is 
A5Vq=  -  AVg-t- AW(- 


where  5pM,  Sp^  where  defined  by  (2-33)  and  (2-34).  Here  we  need  to  point  out  that  the  attraction 
of  a  Bouguer  plate  with  a  variable  density  is  no  longer  a  constant  above  the  Bouguer  plate  and  the 
second  derivatives  of  the  Bouguer  plate  are  no  longer  zero. 

If  the  density  anomaly  5p  is  a  function  of  the  variable  x  and  y,  only,  then  eqs.  (2-33)  and  (2- 
34)  become; 
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8pM  =  h  ■  8p  (x,  V) 

(2-81) 

8p^=|h^Zo  -  jh  j 8p(x,  y| 

(2-82) 

Just  as  in  the  above  derivations  we  give  the  corrections  to  the  first  and  second  derivatives  of  the 
potential  which  were  generated  by  the  density  anomaly  5p; 


5AV,  =  — AVq=- 


q=-jG  Jf  ^apfjdxdy 


8AVy  = 


AVq=4g  f  f  ^Sp^dxdy 

J  J  ^  .  5 


5AV  = 


= - AVq=G  f  f  f  z  Sp  (x,  y,  z)dz dxdy - 

3zp  •'o 

-  — G  f  f  — ^  r  (2zoz-z^)(z(rz)5pdzdxdy 

2  *'  -'t  I  '5  •'0 


Similar  to  eqs.  (2-72)  -  (2-77)  we  have 


8AV.  = 


8AV,y=- 


l^G  ff  5^(x. 

2  L^' 


•'pl(y-yp|  ^•’‘^y 


6AV„  =  . 


=  -3C  f  f  ^  f  iSpdz  dxdy 
-y  •'t  I  '5  *^0 


(2-83) 

(2-84) 

(2-85) 

(2-86) 

(2-87) 

(2-88) 


+  4^  (2z(,z- z^)(zo-z)8pdz  dxdy 


(2-89) 


(2-90) 


r  (zzqZ  -  z^)(z(j-z)  5p  dz  dxdy 
l'7  -^0 


5AV„  = 


=  jGJJ^  f  (- 2zoZ-t- z^jsp  dz  dxdy -t- 

J ~j  (2z(>z- z^)|zo-z)Spdz  dxdy 


(2-91) 


If  the  density  anomaly  is  known  everywhere,  then  we  can  compute  the  change  of  the  potential 
and  its  derivatives  at  the  flight  altitutdc  by  using  the  formulas  as  derived  above.  For  the  airborne 
gradiometry  case,  the  second  derivatives  of  the  potential  of  the  topography  can  be  computed  by 
using  the  Bouguer  reduction  or  the  Heimert's  second  condensation  and  are  added  to  the  measured 
gradient  data  to  obtain  terrain  corrected  gradients.  By  using  the  gradient  data,  the  gravity 
disturbance  (anomalies)  can  be  recovered.  The  disturbing  potential  on  the  ground  then  can  be 
determined  by  using  the  recovered  gravity  data.  The  bias  problem  -  which  occurs  in  the  airborne 
and  satellite  gradiometry,  can  be  solved  by  combining  other  data. 

2.2.3  Indirect  Effect 

Now  we  consider  the  indirect  effect  on  the  potential  of  the  topography  reduction  to  the  geoid. 

Basically,  the  topography  reduction  in  this  repon  is  the  same  as  the  topography  reduction  in 
any  textbook.  Therefore,  the  indirect  effect  is  the  same  as  it  is  described  in  (Heiskanen  and  Moritz, 
19'67,p.  142,  Wichiencharoen,  1982). 

In  the  Bouguer  reduction  the  mass  above  the  sea  level  is  removed  and  the  indirect  effect  is 
very  large. 

The  indirect  effect  of  Helmen's  condensation  is  given  by  (Vanicek  and  Kleusberg,  1987): 

2  2 

2  1  r  ''  h  -hn 

8V,=  -rGph;  -iGp  j  J  -^dxdy 

"  s  (2-92) 


wheres  =  V(x-Xp',^-^(y-y/, 

the  distance  between  computation  point  P  and  current  point  on  the  sea  level. 

Compare  (2-92)  with  (2-49),  we  rind  that  the  change  of  the  potential  due  to  the  Helmen's 
condensation  has  the  same  magnitude  order  as  the  indirect  effect.  If  the  geoid  is  determined 
precisely,  the  indirect  effect  must  be  carefully  considered. 
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2.3  Topographic  Correction  Studied  in  the  Frequency  Domain 

2.3,1  Bouguer  Reduction 

The  derivation  of  the  formulas  for  the  Bouguer  reduction  and  Helmert's  second  condensation 
for  airborne  gradiometry  in  the  space  domain  is  time  consuming.  If  the  reductions  are  studied  in 
the  frequency  domain,  the  problem  becomes  easier. 

We  define  the  2-D  Fourier  transformation  and  its  inverse  as  follows: 


G  (u.v)  =  F  (f'x,  yl) 


If 


f  (X  ,  V 


2ltj(xu  +  yv| 

e  dxdy 


(2-93) 


fix,y|=F  (G(u,v)|=  r  f  G  |u.  vi e ^  *dudv 

(2-94) 

where  F  and  F  '  denote  the  Founer  transformation  and  its  inverse  and  u,  v  are  frequency  variables. 
This  definition  is  a  little  different  from  its  definition  in  most  text  books  (e.g.,  Papoulis,  1968, 
Bracewell,  1965),  but  it  is  more  convenient  for  the  numerical  computation  if  the  subroutines  are 
written  for  the  discrete  Fourier  transformation  which  is  defined  similarly  as  (2-93)  and  (2-94). 

Applying  2-D  Fourier  transformation  to  eq.  (2-1)  we  obtain  the  Fourier  transform  of  the 
potential  Vq  (cf.  Parker,  1972): 


F  {Vo)  =  Ge 


2n:oi 


Y 

n-  1 


iKCi) 


F(p  h 


(2-95) 


here  we  assume  p  =  p  (x,  y),  is  a  function  of  the  variables  x,  y,  and  the  circular  frequency 
a)=(u7-»-v2)l/2;  ZQ  is  the  flight  altitude  of  the  airborne  above  the  sea  level.  If  the  density  p  is 
assumed  to  be  a  constant  it  can  be  taken  out  from  the  Fourier  transformation,  and  eq.  (2-95)  is 
reduced  to  Parker's  solution. 


As  it  was  shown  (Parker,  1972)  that  the  series  convergent  with  a  convergence  rate  better  than 
(hmax/zo)”.  where  hmax  is  the  maximum  nf  the  elevation  of  the  topographs  .  The  physical  meaninc 
of  the  series  is  very  clear,  Wnting  (2-9.31111  the  form 

F  I  Vq)  =  G  — L-  F  { p  (e  ^  I  -  e  j 

2ko) 


=  G 


Ukcd) 


2ko} 


n=  I 


F{p(h-zo)"-p(-  zol")' 


and  taking  the  first  two  terms,  we  get 


(2-96) 
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-F(hp}  +  TiF{(h'-2/,>h|pj 
0)  ' 


(2-97) 


Assuming  that  the  density  of  the  topography  is  constant,  then  eq.  (2-97)  is  nothing  other  than  the 
Fourier  transformation  of  the  potential  of  a  Bouguer  plate  with  thickness  h  plus  a  correction  due  to 
the  variation  of  the  heights.  We  show  this  in  the  following. 

Taking  the  inverse  Fourier  transformation  of  eq.  (2-97)  it  gives 


=  f  f  ^^-dxdy-7tGp(2zoh-h') 

J  J  ^ 

=  Gp  If.  dxdy  -  nGp  (dz^  h-h^j  +  Gp^^  — ~  dxdy 


(2-98) 


where  s  is  the  distance  between  the  computation  point  and  current  point  on  the  sea  level. 

Here  we  have  used  the  relationship: 

and  the  density  p  is  assumed  to  be  a  constant. 

Therefore  the  terms  n>2  in  eq.  (2-96)  represent  the  terrain  effect.  The  terrain  effect  is  smaller 
than  the  effect  of  the  Bouguer  plate,  but  it  is  rough  when  the  topography  is  rugged. 

.Above  computations  are  referenced  to  the  sea  level.  In  order  to  speed  the  convergence  of  the 
senes  (2-95),  reference  surfaces  other  than  sea  level  can  be  chosen.  One  of  the  choices  is  the  mean 
elevation  level,  z  =  hm,  so  that  eq.  (2-95)  becomes 


'  2  n'  \ 


lit  U  i"  •  ■  ,ri  t 

(("■nm)  P  -p  "nil  P  I 


(2-99) 


Becau'.e  we  have 
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^JJ  “  (- hr,  j"S(u,  v|e  ''"*dudv 

=  0,  n  =  1,  2 . 


where  we  have  used  the  property  of  the  5  -  function 


(2-l(X)) 


ff. 


f  (u,  v)  5  (u ,  v)  dudv  =  f  jO,  0) , 


(2-101) 


then  we  can  build  up  an  artificial  formula  to  replace  eq.  (2-98)  if  the  density  p  is  also  a  constant: 
I VJ 

l  Q/c 


F  j  Vgj^  =  Gp  / 1  F  {h !  +  tiF  { h’-2h  h,,,) 


1 


(2 


710) 


Zku) 


n=3 


(ih-M"! 


(2-102) 


The  function  F  (  Vq)c  is  not  exactly  the  Fourier  transformation  of  the  potential  of  the  topography 
with  a  constant  density.  Based  on  (2-l(X)),  the  Fourier  transformation  of  (2-102)  gives  the  same 
potential  Vq  as  the  equation  (2-99)  does.  If  the  density  p  is  not  a  constant,  then  eq.  (2-95)  has  to 
be  used.  The  advantage  of  eq.  (2-102)  is  that  the  series  in  equation  converges  faster  than  the  series 
in  eq.  (2-94). 

The  parameter  h-n  in  (2-102)  is  defined  as  the  mean  elevation  in  the  area.  Based  on  Parker 
(1972)  the  parameter  hm  =  1/2  (h^ax  -  hmin).  where  h^ax  ‘ir'd  ^imin  'S  the  maximum  and  minimum 
elevation  in  the  area,  can  be  the  best  choice  for  the  convergence  of  the  series  in  (2-102).  But,  if  the 
series  converges  very  fast,  the  choice  of  h^  as  a  mean  elevation  or  the  average  value  of  the 
difference  between  the  maximum  and  minimum  elevation  .should  not  play  an  imponant  role. 

In  the  similar  way  as  Parker  did  (ibid)  one  can  show  that  the  convergence  rate  of  the  series  in 
(2-102)  is  better  than  the  .senes 


where  Ff  =  max  ((hrT,.ix  -  hm),  -  (h^in  -  hm)). 

The  vertical  derivative  of  Vn  is  easy  to  get  by  spectral  analysis.  As  in  many  text  boc>ks  (e.g.. 
Bracewell,  1965)  the  relationship  oetween  the  differentia!  operators  and  the  spectra  are  given  in  the 
following: 
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dy 
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=>  JZrtV 


2rtci) 


(2-103) 


k  s  p 

d\  3v  3z 


,  ,  ,p  k^-s/  k  s  P 

(-1)  j  \27c)  U  V  0) 


where  n  =  k  +  s  +  p.  k,  s,  p  =  1,  2,  3....,  and  j  is  unit  imaginary  number. 
The  Fourier  transformation  of  the  attraction  of  the  topography  is  given  by 

F{\M  =  -2rtcuF{VQi 

Insen  eq.  (2-96)  into  eq.  (2-104)  and  taking  the  first  term; 
F(V,}Q=-2rtGF{ph} 

The  inverse  Fourier  transfonnation  of  (2-105)  is 
v;-F'{F{V,)Q}  =  -2;tGph 

It  is  nothing  else  but  the  attraction  of  a  Bouguer  plate  with  thickness  h. 


(2-104) 


(2-105) 


(2-106) 


By  using  eq.  (2-103)  we  get  the  Fourier  transformation  of  the  first  derivatives  of  the  potential 
at  the  flight  altitude: 


FjV,[.j2ttuF;VQj 

F(V,|=j2lvF(V„i 

where  F  (Vq)  is  defined  by  eq.  (?  99)  or  (2-102). 


(2-107) 

(2-108) 


In  the  same  way  we  C)btain  the  Fourier  transformation  of  the  second  derivatives  of  the  potential 
at  the  flight  altitude: 


FiV„)  =  -(2Ku)  F(Vq) 
F(V,^|  =  -(2;t)‘uvF{VQ| 


(2-109) 


(2-110) 
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F|Vy^)  =  -(2TrvrF!VQ| 

(2-111) 

F,'V,,l  =  -j(27r)  uo)F(Vq) 

(2-112) 

F(Vy4  =  -j(27trva)FjVQ) 

(2-113) 

F{V,4=(27ro)fF|VQj 

(2-114) 

Obviously  we  have 

piv  \-i-FIV  '-t-FIV  1=0 

‘  1  ''  xxj  ^  ’  y>|  1  ’  /J.I 

(2-115) 

so  that  there  are  five  independent  second  derivatives  of  the  potential. 

2.3.2  Helmert's  Second  Condensation 

It  is  also  easy  to  apply  the  Fourier  transformation  to  the  Helmert's  second  condensation. 

Applying  the  Fourier  transformation  to  eq.  (2-45),  we  get 

-2na)Z5 

F{Wq}  =  G^ -  F{ph) 

0) 

(2-116) 

This  is  the  Fourier  transfotr.i  of  the  potential  generated  by  the  condensed  topography  layer  on  the 
sea  level.  Obviously,  the  Fourier  transform  of  the  change  of  the  potential  due  to  the  Helmert's 
second  condensation  is  (cf.  eqs.  (2-49)  and  (2-95)): 

f{5v}  =  f{-Vq^Wq) 


=  -  G  - 


1  cL.,"! 


F  \ph 


2:tcu 


If  we  expand  the  exjxment  in  (2-117)  in  series,  then  we  have 
Ir  \  I  ^  ^  (27ro))  (-2ntozo]  f  n\ 

F{5v)  =  -G-i-5l  Z- - ) 

2710) 

The  first  term  of  the  series  in  (2-118)  is  (m^,  n=2) 


(2-117) 


(2-118) 


-G7tF{ph^) 
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and  its  inverse  Fourier  transformation  is 
-  G  7t  p 

which  is  nothing  else  than  the  first  term  in  eq.  (2-49).  Therefore  we  can  expect  that  the  terms 
(n>3)  in  eq.  (2-117)  may  be  ignored  in  most  computations.  This  agrees  with  our  earlier  study  in 
space  domain. 

If  the  topography  is  condensed  on  the  sea  level,  then  the  first  derivatives  of  the  geopotential 
changes  at  the  flight  altitude  by 

5V^=F'‘{j27tuF(5v}} 

SV,  =  F-'{j2mF{8v))  (2., ,51 

5V^=F‘{-27tMF('8v}) 


where  5Vx,  5Vy,  SV^  are  the  change  of  the  first  derivatives  of  the  geopotential,  6V  is  given  by  eq. 
(2-117).  The  effect  of  the  Helmert's  second  condensation  on  the  gravity  gradient  data  is  also  easy 
to  get  by 


5V„=F- 

‘{-(27tufF{5v}/ 

5V,y=F- 

'{-(2n)%v  F  (5v}/ 

5V^,=  F- 

‘{-(2ku)^  f{5v}} 

8V,,=  F 

‘  {-j(27t)^ua)F  (5V}/ 

SV,,=  F- 

'{-j(27t)%0)F{5v}} 

5V,.=  F 

'  (fzTtO))  F  {Sv}) 

In  comparison  with  eqs.  (2-66)  -  (2-71),  eq.  (2-120)  is  simple  and  easy  to  use  for  the  computation. 


3.  Effect  of  Topography  on  Determination  of  the  Geoid  by  Using  Analytical 
Downward  Continuation 

3.1  Introduction 

In  the  preceeding  chapter  we  studied  the  effect  of  the  topography  on  the  airborne  gradiometry; 
some  corrections  are  added  to  the  aerial  gradient  data  to  eliminate  the  topographic  effect.  Now  we 
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consider  the  topography  effect  in  satellite  gradiometry.  This  situation  is  different  from  airborne 
gradiometry  because  a  satellite  flies  much  higher  and  the  high  frequencies  of  the  gravity  field, 
which  mostly  come  from  the  terrain  variation  of  the  earth,  become  small  or  negligable  at  the 
satellite  flight  altitude.  In  this  chapter  we  want  to  find  out  a  topography  correction  which  is  not 
related  to  the  flight  altitude  of  an  airborne  or  a  satellite. 

In  processing  the  satellite  gravity  gradient  data  one  can  ignore  the  mass  above  the  reference 
surface  (e.g.,  the  ellipsoid)  at  first,  the  gravity  disturbance  is  recovered  on  the  reference  surface  by 
using  the  analytical  downward  continuation.  If  such  recovered  gravity  data  are  used  to  determine 
the  coefficients  of  the  spherical  harmonics  of  the  earth's  gravitational  potential,  no  topographic 
correction  is  needed.  This  is  the  same  as  the  analytical  downward  continuation  of  the  free-air 
anomaly  from  the  earth's  surface  to  the  ellipsoid  for  the  determination  of  the  coefficients  of  the 
spherical  hamion'cs  expansion.  If  these  recovered  data  are  used  for  determination  of  the  geoid,  the 
topographic  effect  has  to  be  considered.  The  question  is  what  is  the  difference  between  the 
analytically  downward  continued  potential  and  the  true  one  on  the  ellipsoid.  In  other  words:  What 
is  the  difference  between  the  geoid  determined  by  analytically  downward  continuation  and  the  true 
one? 


The  following  work  is  addressed  on  this  aspect  and  we  focus  our  attention  on  the  topographic 
correction  to  the  geoid  which  is  determined  by  using  the  analytical  downward  continuation  which 
is  the  method  used  in  the  processing  of  the  aerial  gravity  gradient  data. 

3.2.  Remove-Restore  Technique  -  Its  Physical  Meaning  and  Mathematical  Formulation 

In  the  following  we  will  find  the  topographic  correction  to  the  geopotential  which  is 
determined  by  using  the  analytical  downward  continuation  procedure.  In  order  to  do  this,  we 
remove  the  mass  above  the  reference  surface  first,  then  add  it  back.  This  technique  is  called 
remove-restore  technique.  This  technique  is  illustrated  in  the  following. 


Z  Q  Tq  =  Vq  -h  5Tq 


Figure  2.  Geometry  of  computation  of  the  disturbing  potential 


The  disturbing  potential  T  at  the  flight  altitude  Q  consists  of  the  gravitational  potential  (V)  of 
the  mass  above  the  ellipsoid  and  the  potential  (5t)  caused  by  the  density  anomaly  under  the 
ellipsoid. 
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Figure  3.  Removal  of  the  mass  above  the  ellipsoid 

The  mass  above  the  ellipsoid  is  removed.  The  residual  disturbing  potential  6Tq=Tq-Vq  is 
harmonic  above  the  ellipsoid  and  can  be  analytically  downward  continued  to  the  ellipsoid  without 
any  theoretical  problem.  STpQ  is  the  downward  continued  residual  potential  on  the  ellipsoid. 


Z  Q  Flight  altitude 


Figure  4.  Restoration  of  the  mass  above  the  ellipsoid. 

The  mass  above  the  ellipsoid  is  added  back.  The  disturbing  potential  at  point  Pq  located  on  the 
ellipsoid  is  Tp^  =  STp^  +  Vpf,  ,  where  is  the  gravitational  potential  of  the  mass  above  the 
ellipsoid  at  point  Pq. 
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Now  we  assume  that  the  disturbing  potential  Tq  is  analytically  downward  continued  to  the 
ellipsoid;  Tp„  =  5TpQ  +  V* 


*Pq  —  V/ipQ  ■  TpQ^Viril^lW  J 

Vq  at  the  flight  altitude,  respectively. 


where  Tjfg  and  Vp^^  are  the  downward  continued  potentials  of  the  Tq. 


r  * 

'Po 


Compare  Tp  with  the  disturbing  potential  Tp,  which  is  the  original  disturbing  potential  on  the 
ellipsoid,  we  get 


8V=Tp„-rp,  =  8Tp„  +  Vp„-8Tp,-v;„ 


=  -  V  +  V 

Po  Po  (3.1) 

Obviously,  8V  defined  by  (3-1)  is  the  difference  between  the  true  gravitational  potential  of  the 
topographic  mass  at  point  Pq  and  an  anificial  potential  Vp^  which  is  the  gravitational  potential  of 
the  topographic  mass  at  point  Q  (outside  the  Earth’s  surface)  analytically  downward  continued  to 
point  Pq. 

From  the  above  short  description  of  the  remove-restore  technique  we  can  see,  8V  is 
independent  of  the  flight  altitude;  therefore  it  is  suitable  for  computation  of  the  topographic  effect 
on  the  determination  of  the  geoid  by  using  analytical  downward  continuation.  For  instance,  the 
disturbing  potential  Tq  can  be  computed  by  processing  the  satellite  or  the  airborne  gravity  gradient 
data  analytically  downward  continued  onto  the  ellipsoid.  So  the  disturbing  potential  Tp^  can  be 
obtained  by  adding  the  5V,and  at  the  same  time  the  geoid  is  also  determined  from  Brun's  equation. 

We  still  call  SV  the  topo^aphic  correction  because  it  comes  from  the  earth's  topography.  If 
the  eanh  is  exactly  an  ellipsoid  which  coincides  with  the  reference  ellipsoid,  then  8V  is  equal  to 
zero. 


Before  we  compute  the  topographic  correction  8V,  we  first  consider  the  correction  of  a 
spherical  Bouguer  plate. 


3.3  Effect  of  A  Spherical  Bouguer  Plate  on  the  Analytical  Downward  Continuation 

We  approximate  the  ellipsoid  by  a  sphere  with  the  radius  R,  the  mean  radius  of  the  earth.  A 
spherical  Bouguer  plate  has  the  boundanes  r  =  R  and  r  =  R  -t-  hp,  hp  is  the  thickness  of  the 
Bouguer  plate.  ^ 

The  geopotential  of  the  Bouguer  plate,  just  like  a  point  mass  or  a  homogeneous  sphere,  is 
given  by 


V 


Q 


CM 

R+Zo 


(3-2) 


where  Vq  is  the  geopotential  of  the  Bouguer  plate  at  the  flight  altitude  and  zq  is  the  flight  altitude 
above  the  ellipsoid;  G  is  the  gravitational  constant  and  M  is  the  mass  of  the  Bouguer  plate  which  is 
given  by 


M  =  i7rp[(R  +  h/ 


-  R 
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=  47tp 


R^p+RhJ  +  ihJ 


(3-3) 


where  we  have  assumed  that  the  density  of  the  Bouguer  plate  is  a  constant. 

If  we  analytically  downward  continue  the  potential  Vq  from  the  flight  altitude  to  the  ellipsoid, 
we  obtain 


GM 

R 


=  47tGp 


(3-4) 


where  is  the  downward  continued  potential  on  the  ellipsoid  at  Pq. 

At  the  same  time,  the  gravitational  potential  of  the  Bouguer  plate  at  point  Pq  is  given  by 


Vp,=  2:tGp 


|R  +  h. 


=  47tGpR  hp+ 2;:Gphp  ^3  5^ 

The  mean  radius  of  the  eanh,  R,  is  much  larger  than  the  thickness  of  the  Bouguer  plate.  For 
example,  if  we  take  the  maximum  thickness  of  the  Bouguer  plate  as  hmax  =  10  km,  then  the  ratio  of 
the  thickness  of  the  Bouguer  plate  and  the  mean  radius  of  the  earth  is 


10 

6.4  X  10^ 


=  1.6  X  10'^ 


(3-6) 


Therefore  the  third  term  (ZaGphg)  in  (3-5)  is  one  thousand  times  smaller  than  the  second  one. 
Comparing  (3-4)  with  (3-5),  we  find  the  potential  difference 


5V  =  - 


V=-2^Gph 


p 


(3-7) 


Trie  term  on  the  order  of  hp/R  in  5V  is  neglected. 

Therefore  the  analytically  downward  continued  geopotential  of  a  Bouguer  plate  from  outside 
onto  the  ellipsoid  differs  from  the  original  geopotential  by  -  27tGp  h^  This  result  has  been  known 
from  many  previous  publications  (e.g.,  Wichiencharoen,  1982,  p.  16,  Table  1).  Here  we 
emphasize  the  fact  that  the  potential  difference  -2hGphp  is  due  to  the  use  of  the  analytical 
downward  continuation. 

3.4  Topographic  Effect  on  the  Determination  of  to  the  Disturbing  Potential  on  the  Ellipsoid  by 

Using  Analytical  Downward  Continuation 

Now  we  consider  the  topographic  effect  on  the  disturbing  potential  which  is  analytically 
downward  continued  from  outside  earth's  surface  onto  the  ellipsoid.  The  geometry  of  the 
computation  of  the  potential  of  the  topography  is  drawn  in  Figure  5. 
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Figure  5.  Geometry  of  the  Computation  of  the  Potential  of  the  Topography 

The  gravitational  potential  of  the  mass  above  the  ellipsoid  at  the  flight  altitude  is  given  by  (cf. 
Moritz,  1968); 


where  the  density  p  is  assumed  to  be  a  constant,  L(p  is  the  distance  between  the  point  Q  and  the 
current  point  inside  the  topography;  O  is  the  unit  sphere  and  R  is  the  mean  earth  radius. 

The  potential  Vq  can  be  separated  into 

Vq=V, +V2  (3-9) 

with 


V,  =  GpR'rr 

*'o  •'o  ^<0 

V2  =  GpR'rr  /'^do 


(3-10) 

(3-11) 


where  V ]  is  the  potential  of  the  spherical  Bouguer  plate  and  V2  is  the  effect  of  the  terrain  around 
point  P.  Without  loss  of  accuracy  (3-11)  can  be  approximated  by  a  plane  approximation 


(3-12) 
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where  T  is  the  tangential  plane  of  the  ellipsoid  at  point  Po.  and 

L- =  S“  +  (ZQ  -  z)“  (3-13) 

where  L  and  s  were  previously  defined  (p.  3).  Now  we  analytically  downward  continue  Vq  to  the 
ellipsoid  to  get  the  potential  difference  5V.  Based  on  last  section  the  potential  difference  of  a 
Bouguer  plate  is  -27tGphp  In  the  following  we  consider  only  the  terrain  effect  V2. 

We  assume  that  the  potential  V2  is  analytically  downward  continued  to  the  ellipsoid  by  the 
Taylor's  series: 


V2=V2-Zo 


avj 

dZy 


2 

0 


2 

d  Vj 


1  3  3  V2 

6  3 

3z„ 


For  this  series  we  have 


(3-14) 


a  1  Zq-z 


a‘  1 


a  1 

3  L 
aZy 


(3-15) 


(3-16) 


(3-17) 


where  0  (L‘^)  denotes  the  term  contains  L‘^.  Putting  (3-15),  (3-16),  (3-17)  and  (3-12)  into  (3-14) 
we  get 


v:  =  Gpff  r'^dxdy 

•'r  hp  L  (3-18) 

with 

2  -  3  2  2 

_L-1  ^0'  ^^0^  3  -  Zpz-t-  ZyZ 

,  *  “  L  ^  01  3  ^2  ,5 

I  2L  L  (3-19) 


Because  we  have  always 
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S  +  Zn  >  2Zn  z  -  z 


or 


s  +  z, 


one  can  write 


1  2zoZ  - 

2  ,  '2 


2 


(3-20) 


(3-21) 


I  3  '3 

L  L 


1  + 


3  2zoZ  -  z‘ 


with 


,'2  2^2 

L  =  s  -t-  Zq 

Insert  (3-21),  (3-22)  and  (3-23)  into  (3- 19),  we  obtain 

2  2  ,2  4' 

1  1  Z(,-Z  3  -2z,>7.  +  z 

_  i  + - +  _ - 


(3-22) 

(3-23) 


(3-24) 


(3-25) 


and  the  potential  V7  is  given  by 
■  h 


i,  *  jGp// 


dzdxdv  + 


■f’  -2z,V.^+z'' 


dzdxdy 


(3-26) 


Here  we  need  to  point  out  that  the  quantity  V2  has  no  physical  meaning.  If  we  say  it  has  a 
meaning,  it  means  that  it  is  a  fictitious  potential  created  by  analytically  downward  continuing  the 
potential  of  the  terrain  from  outside  the  eanh  onto  the  ellipsoid. 

The  potential  of  the  topography  at  point  Pp  is  given  bv 
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where  is  the  distance  between  the  point  Pp  and  the  current  point  ins.Je  the  topography. 
The  potential  Vp^  can  be  split  into 


(3-27) 


V  =  V  -t-  V 

P'  Po  Po 


(3-28) 


with 


v,,=  Gpr7/^  X 


'p  dz 


VX^=GpR^//^  —  da 


(3-29) 


(3-30) 


The  potential  Vp^  is  the  potential  of  the  spherical  Bougucr  plate,  but  the  computation  point  po  is 
under  the  Bouguer  plate  on  the  ellipsoid.  Therefore  (3-29)  is  the  same  as  eq.  (3-5).  For  the 
potential  Vp^  we  take  the  plane  approximation; 
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-^/xx: 


dz 
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with 


c  s  ,  "I 
b  -  =  S'  -1-  z^- 


(3-31) 


(3-32) 


We  expand  (3-31)  into  a  senes: 

^7/  X ,  X 7‘7X  X 


'‘’ff.i.T. 


■  dzdxdy 


putting  (3-33)  and  (3-26)  together,  we  get 


(3.331 
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(3-34) 


Completely  we  have  the  potential  difference  5V  is  a  sum  of  the  5V2  and  a  term  -ZrtGphp  which  is 
from  the  downward  continuation  of  the  potential  through  a  Bouguer  plate: 


5V  =  -27tGphp 


h-h. 


-  dxdy 


(3-35) 


The  potential  difference  5V  represented  by  (3-35)  is  a  function  of  the  flight  altitude  zo.  This  is  in 
contradiction  with  our  earlier  statement  that  the  potential  difference  5V  is  independent  of  flight 
altitude  zq.  This  contradiction  comes  from  the  approximation  used  in  the  derivation  of  5V.  In 
frequence  domain  the  problem  can  be  studied  without  approximation,  and  we  can  get  the  5^’ 
independent  of  flight  altitude  as  expected  (cf.  next  section). 

Equation  (3-35)  gives  the  potential  difference  in  space  domain  approximately.  For  this 
equation  we  have  two  more  things  to  say: 

1.  The  expression  of  5V  is  an  approximation,  the  integrals  which  contain  terms  higher  than  S'5 
were  omitted.  The  flight  altitude  is  assumed  above  the  topography.  The  above  derivations  are  also 
valid,  if  the  computation  point  is  chosen  on  the  topography  of  the  Earth.  In  this  case,  eq.  (3-35) 
becomes: 


6V  =  .2llOphs  *  1  Cphl  (  I  dxdy 

*  if  ,,  ,d) 

2.  The  last  term  in  (3-36)  is  a  very  small  quantity,  A  numerical  test  (Wang,  1989)  showed  that  the 
last  term  in  (3-36)  contributed  to  the  geoid  at  the  millimeter  level  in  a  rough  mountain  area  so  that  it 
can  be  neglected  in  most  ca,se.s; 


The  change  of  the  geoid  caused  by  5V  is  then  given  by  the  Brun  s  fonnula: 


Yo 


(3-37) 


where  yo  is  the  normal  gravity  on  the  ellipsoid 

The  first  term  in  (3-36)  is  primary  and  its  contribution  to  the  geoid  can  reach  1-2  meters  in 
high  mountains,  therefore  the  topogra'phic  correction  is  significant  and  cannot  be  neglected. 

Now  we  consider  the  effect  of  the  topography  to  the  deflections  of  the  venical. 


It  is  well  known  that  a  Bouguer  plate  with  constant  density  has  no  contribution  to  the 
deflection  of  the  vertical.  Here  an  implicit  assumption  is  made  that  the  elevation  hp  is  a  constant. 


34 


Therefore  hp  c.tn  be  viewed  as  a  constant  in  eq.  (3-36)  too,  when  the  denvatives  are  applied,  Th 
corrections  to  the  detlections  of  the  vertical  due  to  the  5V  are  given  then 


1  36V 

Yo 


Yo 


1 


( X  -  Xp)  dx  dy 


_  1  36V 

5n  - - 

Yo  3yp 


(3-3S) 


i>  ->pi<ixdy 


(3-39) 


3.5  Topographic  Correction  in  the  Frequency  Domain 

Now  we  study  the  topn-igraphic  correction  studied  in  frequency  domain. 

The  geometry  of  the  computation  of  the  topographic  correction  is  drawn  in  Figure  6. 


2 


A 


Figure  6  Geometry  of  the  computation  of  Topographic  Correction 
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In  the  plane  approximation  the  gravitatior  >1  potential  ot  the  topography  at  points  Q  and  Po  is 
izivcn  bv: 


,=  G/r  pl..y|/'’!|idxjy 


(3-40) 


r°fj  Plx.yl  r'ldxdy 


(3-41) 


Here  we  have  assumed  the  density  of  the  topography  is  a  function  of  the  the  horizontal  variables  x, 
y,  and  the  distances  between  the  current  point  and  the  computation  points  Q  and  Pp  are 

L-  =  s‘  -t-  (ZQ  -  z)- 

S  ^  =  s-  +  z- 

They  were  previously  defined. 

Applying  2  D  Fourier  transformation  to  (3-40),  we  get  (cf.  Schwiur.  et  al.,  1989): 


F{Vq}-G - ;F(p(e  '^'-e  If 


=  G  e  - -  F  jp  le  ■  ll, 

2:t(o 


(3-42) 


Analytical  downward  continuing  the  jxiteniial  Vgonto  the  ellipsoid,  then  its  Fourier  transformation 
is  given  by: 


(3-43) 


PiVQ|=G - ^  F^p(e 


In  the  same  manner  we  obtain  the  Fourier  transfonnaiion  of  Vj 
„  .  .  ,  ^  1  I  I  2ni:)li  1  ^ 

!^Po/  =  -  -  -y  f-  'fU’  -  ii; 


(3-44) 


Combine  eqs,  (3-43 1  and  (3-44)  we  get  the  F-’ouner  transform  of  the  potential  difference  5V: 


f(6v;  =  f|- VpJ 
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=  -  G 


1  p  / 

- 
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(3-45) 


Eq.  (3-45^  is  independent  of  the  altitude  zo  as  expected.  Expanding  (3-45)  into  a  series,  then  we 
get 


F  {5v}  =  -G 


TTO) 
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(2110)) 


F\ph  j 


.  2n 


=  -2GFi.,h')-G-4f 


710) 


The  inverse  Fourier  transformatioi  of  the  first  term  of  eq.  (3-46) 
-  2kG  F  {ph2) 


(3-46) 


(3-47) 


is  equal  to  -ZTCjph^  which  is  nothing  other  than  the  correction  of  a  Bouguer  plate  (cf.  eq.  (3-7)). 

Clearly  the  magnitude  of  the  series  in  (3-46)  depends  very  strongly  on  the  roughness  of  the 
lOpography:  if  the  topography  is  smooth,  the  series  has  little  contribution  to  the  change  of  the 
potential  oV,  if  h  =  constant,  the  senes  in  (3-46)  has  no  contribution  to  5V. 

The  last  sum  in  (3-46)  supplies  mostly  very  high  frequencies  in  5V,  therefore  it  is  not 
important  to  6V  which  is  assumed  that  the  very  high  frequencies  have  few  contribution.  This 
agrees  with  our  earlier  study  in  space  domain. 

If  we  compute  (3-46)  and  its  inverse  Fourier  transformation,  the  discrete  Fourier 
transformation  has  to  be  used.  In  this  ca.se  the  influence  of  aliasing  and  leakage  becomes  a  very 
serious  problem,  because  both  phenomenan  have  the  most  effect  on  high  frequencies.  If  such 
wrong  high  frequencies  are  amplified  by  O)”,  then  results  will  not  be  reliable.  Therefore  the  low- 
pass  filters  is  needed  for  such  computation. 

The  first  derivatives  of  6V  is  given  by 

5V,  =  j27tF'’  {uF{5v)} 

6V^  =  j27iF  '  {vf{5v)} 

5V^=-27tF  ‘(  wF  {5V}} 


(3-47) 
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where  F  (5V)  is  given  by  eq.  (3-46). 

In  the  above  discussion  we  have  not  considered  whether  the  analytical  downward  continuation 
is  convergent.  We  believe  the  above  derivation  is  justified  based  on  the  assumption  that  the  series 
in  eqs.  (3-46)  is  at  least  numerical  convergent.  We  cannot  guarantee  that  the  series  in  eq.  (3-45)  is 
convergent.  The  terms  (n>2)  in  the  series  may  be  very  rough,  but  they  could  be  very  small  after 
some  smoothing.  Ignoring  the  terms  (n>2)  in  eq.  (3-45)  may  not  cause  serious  error  in  the 
computation. 

In  the  next  section  we  will  study  the  convergence  problem  of  the  analytical  downward 
continuation  under  planar  approximation. 

3.6  Convergence  Study  of  the  Analytical  Downward  Continuation  Under  Planar  Approximation 

Even  though  the  convergence  problem  of  analytical  downward  continuing  a  harmonic  function 
into  the  mass  does  not  have  too  much  meaning  on  the  practical  point  of  view,  it  is  still  interesting  to 
know  the  mathematical  property  of  the  analytical  downward  continuation  under  planar 
approximation. 

Equation  (3-46)  is  directly  related  to  the  convergence  or  divergence  of  the  analytical 
downward  continuation  of  the  geopotential  inside  the  earth.  If  (3-46)  is  convergent,  then  the 
analytical  downward  continuation  is  also  convergent.  Obviously  the  convergence  problem 
depends  very  strongly  on  the  roughness  of  the  topography.  If  h=0,  then  F  (6V)=  0;  if  h  = 
constant,  then  F"'  fF{5V})  =  -  27tGph‘,  the  topographic  correction  of  a  Bougucr  plate  with 
thickness  h.  In  both  cases  we  do  not  have  any  divergence  problem. 

From  eq.  (3-45)  we  have 


F^ 

'  ”  2 (n+l)(n+2)...2n 

TtO) 

n!  / 

(3-48) 


we  assume  that  the  Fourier  transformation  of  ph'F/n'  exists,  it  means,  the  following  condition 
holds: 


2n 


j  I  ^ =  M.  <  - 

From  (3-49)  one  can  get 

2n+2  2  2n 


n-t-l 


(3-49) 


(3-50) 
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Thir.  means  that  the  following  inequality  holds; 
Mn+1  hL, 


Mn  n+1 


For  the  Fourier  transformation  of  the  quantity  ph-^/n!  we  have 

In 


(3-51) 
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2nj  (xu+yv|  ph 


n! 


dxdy 
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(3-52) 


therefore  for  eq.  (3-48)  we  have 


|F{5v}|<-^y _ _ 
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(3-53) 


The  convergence  radius  of  the  series  in  (3-53)  is  easy  to  get  by  using  D'Alembert's  criterion 
{  n  +  1)  (n  +  2)...2n 
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(3-54) 


Therefore  the  Fourier  transformation  of  the  potential  difference  potential  6V  is  convergent  in  the 
frequency  domain  0<aK«>.  Note  that  the  infinitive  point  ci>=»  is  not  included  in  the  convergence 
domain.  In  the  same  way  we  can  show  that  the  inverse  Fourier  transformation  of  eq.  (3-46)  is 
convergent,  provided  'he  infinite  point  a)=«>  is  not  included  in  F  (SV). 

We  assume  the  maximum  frequency  of  the  Fourier  transformation  of  6V  is  cuo,  where  coo  can 
be  an  arbitrarily  large  number,  but  not  equal  to  infinity.  From  eq.  (3-46)  we  have 
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(3-55) 


Here  we  have  used  the  polar  coorcinates  system  in  the  frequence  domain.  The  relationship 
between  the  frequences  variable  u,  v  and  the  polar  coordinates  system  is  defined  as: 


(  u  =  CO  cos  a 
I  V  =  CO  sin  a 


The  convergence  radius  of  the  series  in  (3-61 )  is  given  by 


d  =  lim 
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(3-55a) 


(3-56) 


From  eq.  (3-56)  we  know  that  the  series  in  (3-55)  is  convergent  for  all  coq:  0<coo<«>. 

From  the  above  derivation  we  can  say  that  the  Fourier  transformation  of  the  6V  and  the 
potential  difference  5V  are  convergent  everywhere,  except  at  the  point  co='».  That  means  the 
analytical  downward  continuation  (ADC)  of  a  harmonic  function  into  the  earth  is  also  convergent, 
except  at  the  point  co=<».  This  result  is  slightly  a  surprise,  because  the  ADC  may  be  assumed 
(considered)  as  divergent  (Sjoberg,  1977).  We  have  shown  that  the  analytical  downward 
continuation  of  the  external  potential  into  the  earth  is  almost  convergent  everywhere!  Here  "almost 
everywhere"  means  the  infinite  point  Cj»=oo  is  not  included. 

Indeed,  the  convergence  of  analytical  downward  continuing  the  geopotential  outside  the 
earth's  surface  into  the  mass  of  the  earth  is  a  delicate  problem.  From  eqs.  (3-43)  or  (3-46)  one  can 
see  that  the  convergence  problem  depends  on  the  roughness  of  the  topography.  Based  on  the 
above  derivation  we  know  that  the  series  in  (3-46)  is  covnergent  for  all  to:  At  the  infinite 
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point  co=«>,  the  convergence  of  (3-46)  is  not  defined.  Obviously,  the  series  is  convergent  at  the 
infinite  point  only  if 


lim 


27t(o)""’^F 


<oo 


(3-57) 


exists. 

Note  that  the  Laplace  operator  in  spectral  domain  is  corresponding  to 

3^  ,2 

A:= - - —  (2710)) 

ax^  V 

the  eq.  (3-57)  is  equivalent  to 

lim  F{(-l)"A'’‘'(ph^")}<oo 


(3-58) 


(3-59) 


where  A"'*  =  AAA  A  (n  -  1  times).- 

Eq.  (3-59)  indicates  that  the  function  ph^"  must  be  an  analytical  function  (an  analytical  function  has 
infinite  derivatives).  This  is  not  the  case  for  the  real  earth,  because  there  are  steep  cliffs  and 
overhanging  rocks  on  the  earth's  surface.  In  this  meaning  the  ADC  is  devergent. 

But  one  has  to  keep  in  mind  that  the  real  world  can  only  be  approximated  by  some 
mathematical  models.  If  the  elevation  h  is  approximated  by  an  analytical  surface,  or,  more 
practical,  be  approximated  by  point  or  mean  blocks  values,  whose  Fourier  transformation  has  finite 
frequencies,  then  (3-46)  is  convergent.  Let  us  now  imagine  in  the  future;  during  the  developement 
of  the  technology  the  topography  of  the  earth  is  measured  more  and  more  in  detail,  we  get  more 
and  more  elevation  data  which  includes  high  and  very  high  frequencies  of  the  topography.  But  the 
number  of  the  data  is  always  limited  and  its  Fourier  transformation  has  always  finite  frequencies. 
In  this  case,  eq.  (3-46)  is  always  convergent.  Therefore  we  will  sa/  that  the  ADC  is  convergent 
for  all  practical  use. 

The  above  conclusion  is  also  imponant  for  gcopotential  modeling.  If  the  ADC  is  convergent, 
we  can  expand  the  spherical  harmonics  of  the  gcopotential  up  to  high  and  very  high  degree  and 
order  without  any  theoretical  problem.  In  practice  the  geopotential  can  be  modeled  with  a  finite 
degree  and  order  which  correspond  to  a  certain  frequency  to  =  (Oq,  therefore  the  ADC  can  be 
considered  as  convergent  for  modeling  the  spherical  harmonic  of  the  geopotential. 
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4.  Effect  of  the  Topography  on  Solving  Molodensky's  Problem  and  the 
Determination  of  the  Coefficients  of  the  Spherical  Harmonics  of  the 
Geopotential 

4.i  L'.roducuon 

In  solving  Molodensky's  problem  and  the  determination  of  the  coefficients  of  the  spherical 
harmonics  of  the  gravitational  potential  of  the  earth  the  analytical  downward  continuation  (ADC)  is 
widely  used.  Moritz  (1980,  section  45)  used  ADC  to  continue  the  gravity  anomaly  to  the  point 
level  and  the  disturbing  potential  is  determined  by  using  the  Stokes’  integral.  Rapp  et  al.,  (1986) 
also  used  ADC  to  continue  the  gravity  anomaly  to  the  ellipsoid  for  the  determination  of  the 
coefficients  of  the  spherical  harmonics.  It  is  always  assumed  that  ADC  is  correct  for  the 
determination  of  the  disturbing  potential  on  the  earth's  surface.  The  disturbing  potential  on  the 
eanh's  surface  and  outside  the  earth  can  also  be  represented  by  the  spherical  harmonic  expansion 
determined  by  using  ADC  of  the  gravity  anomaly  to  the  ellipsoid.  In  this  chapter  we  will  look  at 
this  problem  in  more  detail.  We  wish  to  answer  the  questions:  Can  the  spherical  harmonic 
represent  the  disturbing  potential  on  and  outside  the  earth  exactly?  This  problem  is  different  from 
the  divergence  or  convergence  of  the  spherical  harmonics  on  the  eanh's  surface.  Based  on  the 
derivation  in  chapter  3  we  assume  that  an  spherical  harmonic  expansion  of  the  geopotential,  such 
as  used  by  Rapp  (1984,  eq.  (1)),  be  convergent  on  the  eanh's  surface  and  even  so  on  the  ellipsoid. 
This  assumption  may  be  argued.  But  we  think  it  is  reasonable  for  the  practical  use.  The 
convergence  or  divergence  of  the  series  is  meaningless  in  practice  because  a  small  "grain  of  sand" 
can  make  a  series  divergent  (Moritz,  1980,  p.  64),  and  the  potential  of  a  grain  of  sand  is  so  small 
that  it  can  never  be  included  in  a  gravity  model  related  to  the  eanh's  gravity  field.  Now  the 
question  is:  does  the  expansion  (ibid)  represent  the  gravity  field  of  the  eanh  on  and  outside  the 
eanh's  surface  well?  This  problem  is  imponant  for  today's  gravitational  modeling.  Up  to  now  we 
have  various  geopotential  models  (e.g.,  GEMTl,  OSU89. ...)  generated  by  different  groups  and 
by  using  different  data  sets.  If  the  spherical  harmonics  defined  by  (Rapp,  1984,  eq,  (1)) 
represents  the  gcopotential  on  and  outside  the  eanh's  surface  exactly,  then  the  geopotential 
modeling  can  be  expanded  into  high  and  very  high  degree  and  order.  An  explicit  example:  if  the 
spherical  harmonics  is  expanded  in  degree  and  order  10000,  then  the  geopotential  is  modeled  with 
resolution  2  km.  This  model  will  represent  a  global  uniform  detailed  gravity  field.  But  if  the 
spherical  harmonics  (ibid)  do  not  represent  the  geopotential  on  and  outside  the  eanh's  surface  well, 
then  such  expansion  is  meaningless,  because  the  high  degree  and  order  of  the  spherical  harmonics 
may  be  totally  wrong. 

In  practice  the  high  degree  and  order  of  the  spherical  harmonics  are  effected  mostly  by 
different  errors  and  are  not  as  reliable  as  the  lower  degree  and  order.  But  if  the  expansion  (ibid) 
gives  the  geopoteniial  on  and  outside  the  eanh’s  surface  as  well  as  we  wish,  one  can  expect  that  the 
more  and  more  gravity  data  will  be  gathered  with  more  and  more  better  accuracy  and  data 
coverage,  during  the  development  and  improvement  of  the  high  technology  and  space  technology, 
and  an  expansion  of  spherical  harmonics  to  high  degree  and  order  will  approximate  the 
geopotential  on  and  outside  the  eanh's  surface  to  a  great  extent. 

In  the  following  we  will  study  the  effect  of  the  topography  on  solving  the  Molodensky's 
problem  and  determination  of  the  coefficients  of  the  spherical  harmonics  of  the  earth's  gravitational 
potential.  More  exactly,  to  study  the  error  of  the  analytical  downward  continuation  for  solving 
Molodensky's  problem  and  the  determination  of  the  coefficients  of  the  spherical  harmonics.  A 
similar  study  on  "error  of  the  analytical  downward  continuations"  was  given  by  Sjoberg  (1977, 
1980).  The  error  was  studied  on  a  sphere  earth  and  the  spherical  harmonics  were  used  for  the 
study.  Here  we  will  study  this  problem  for  a  flat  earth.  The  planar  approximation  should  not 
cause  serious  problem  for  our  study.  Note  that  we  use  the  term  "effect  of  the  topography"  instead 
the  "error  of  the  downward  continuation",  in  the  following. 
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4.2  Analytical  Downward  Continuation  for  Solving  Molodensky's  Problem 

Anc  was  used  for  colving  the  Molodensky's  problem  (Moritz,  1980  sertion  45).  This 
method  corresponds  to  the  analytical  continuation  of  the  externa!  potential  f  into  the  earth's  interior 
(ibid,  p.  378)  on  the  point  level. 

Now  let  us  look  at  this  problem  in  more  detail.  This  method  is  illustrated  in  Figure  6. 


Figure  7.  Geometry  of  Solving  Molodensky's  Problem  by  Using  Analytical  Downward 
Continuation 

Gravity  data  are  given  on  the  earth's  surface  and  then  analytically  downward  to  the  point  level. 
The  disturbing  potential  at  the  point  P  is  then  determined  by  using  Stokes'  integral.  This  method  is 
equivalent  to  the  analytical  downward  continuation  of  the  external  potential  down  to  the  point  level 
so  that  the  effect  of  the  topography  in  this  procedure  can  be  obtained  from  the  derivation  in  Chapter 
3. 


Based  on  the  description  in  section  3.2,  the  difference  between  the  disturbing  potential 
determined  by  using  the  ADC  and  the  original  one  is  given  by 


(4-1) 


where  Vp  are  the  potential  determined  by  using  ADC  and  the  potential  of  the  mass  above  the 
point  level  respectively. 

Comparing  Figure  7  with  Figure  6,  we  get  the  Vp  similar  to  eqs.  (3-43)  and  (3-44); 
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Combining  (4-2)  and  (4-3)  we  obtain 
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(4-5)  shows  the  difference  between  the  potential  obtained  by  using  ADC  and  the  true  one. 
Obviously,  (4-5)  is  the  effect  of  the  topography  which  is  neglected  in  the  solving  of  the 
Molodensky's  problem  by  using  the  ADC. 


Theoretically,  (4-5)  plays  a  role  of  the  convergence  of  the  ADC  of  a  harmonic  function  onto 
the  eanh.  It  is  clear  that  the  convergence  of  the  scries  in  (4-5)  depends  on  the  roughness  of  the 
topo^phy.  As  is  shown  in  section  3.5,  one  can  show  that  the  topographic  correction  5Vp  and  its 
Fourier  transformation  are  convergent  almost  everywhere  except  at  the  infinite  point  of  the  circle 
frequency  co  =  <».  Obviously,  the  5Vp  could  be  rough,  but  if  some  kinds  of  smoothing  are  taken. 
(4-5)  could  be  very  small  and  be  neglected.  This  is  also  the  reason  why  we  assume  that  the  ADC 
is  always  correct  for  the  practical  using.  Approximately  we  take  the  first  two  terms  in  (4-5): 
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For  point  P  on  the  topographic  surtace.  the  first  term  has  no  contribution  to  6Vp.  We  have  then 
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This  term  should  be  very  small  after  the  smoothing  and  can  be  omitted. 


(4-7) 


Therefore  we  can  say  The  using  t>f  the  ADC  to  solve  the  Molodensky's  problem  is 
reasonable  and  correct.  The  ti^pography  effect  6Vp  in  this  procedure  represents  very  high 
frequency  in  the  disturbing  potential  and  m  practice  it  can  be  neglected. 
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4.3  Analytical  Downward  Continuation  for  Determination  of  the  Spherical  Harmonics  of  the 

Earth's  Gravitational  Potential 

If  wc  aiialy  i.ic<iHy  downward  continue  tne  gravity  anomaly  to  the  ellipsoid  and  this  reduced 
gravity  anomaly  is  used  for  the  determination  of  a  spherical  harmonic  function,  we  always  assume 
that  this  spherical  harmonic  function  represents  the  disturbing  potential  on  and  outside  the  earth 
very  well.  In  the  following  we  will  show  that  this  spherical  harmonic  function  gives  the  disturbing 
potential  above  and  on  the  Briilotiin  sphere  (Moritz,  1980,  p.  431)  exactly  if  this  series  is 
convergent.  There  is  an  effect  of  the  topography  in  the  vicinity  of  the  eanh's  surface  (the  space 
between  the  topography  surface  and  the  Brillouin  sphere)  to  the  geopotential  represented  by  this 
spherical  harmonic  expansion.  This  means  that  a  single  splierical  harmonic  expression  for  the 
representation  of  the  disturbing  potential  outside  the  earth's  surface  is  not  enough.  The  attempt  to 
include  the  topographic  effect  in  a  series  for  the  representation  of  the  disturbing  potential  in  the 
whole  space  can  be  found  in  (Petrovskaya,  1976,  1977).  This  method  seems  not  so  promising  for 
the  practical  use. 

Here  we  assume  we  have  a  sphencal  harmonic  expansion  determined  by  using  ADC  and  we 
intend  to  find  the  topographic  effect  near  by  the  earth's  surface. 

In  the  first  place  we  will  show  that  if  the  analytical  downward  continuation  of  the  spherical 
harmonics  to  the  ellipsoid  exists  (convergent),  and  if  the  point  P  is  on  the  Brillouin  sphere,  then  the 
downward  continued  spherical  harmonics  represent  the  disturbing  potential  at  point  P  exactly. 

We  denote  the  radius  of  the  Brillouin  sphere  by  Rb,  and  the  semimajor  axis  of  the  ellipsoid  by 
a.  On  the  Brillouin  sphere  the  spherical  harmonic  expansion  of  the  disturbing  potential  is 
convergent  (no  mass  above  the  Brillouin  sphere): 


^  X  Z  (Cnm  cos  mX  +  s  sin  mx)  Pnn,('sin  (p) 
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where  Tg  is  the  disturbing  potential  on  the  Brillouin  sphere,  kM  is  the  gravitational  constant  times 
the  earth  s  mass;  Cnm  and  Snm  are  fully  normalized  potential  coefficients.  If  we  analytically 
downward  continue  the  Tg  to  a  sphere  with  radius  a,  then  we  have 
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It  is  difficult  to  say  whether  the  series  in  (4-9)  is  divergent  or  covergent.  If  the  series  in  (4-9)  is 
convergent,  on  the  Brillouin  sphere  we  have  r=RB.  and  (4-9)  returns  to  (4-8). 

For  the  points  above  the  Brillouin  sphere  we  have  the  same  conclusion. 

Therefore,  the  spherical  hanmonic  expansion  determined  by  using  the  ADC  method  represents 
the  disturbing  potential  on  and  above  the  Brillouin  sphere  (the  bounding  sphere)  exactly,  provided 
the  downward  continued  potential  is  convergent.  Even  if  the  series  in  (4-9)  is  divergent,  based  on 
the  Runge's  theorem,  we  can  find  a  spherical  harmonic  expansion  similar  to  (4-9)  to  approximate 
the  disturbing  potential  as  well  as  we  wish  (cf.  Moritz,  1980,  p.  67). 

Now  we  are  going  to  find  out  the  difference  betwe»n  the  true  gravitational  potential  at  point  P 
and  the  gravitational  potential  represented  by  a  spherical  harmonic  expansion  like  eq.  (4-9).  If 
point  P  is  on  the  earth's  surface,  the  point  level  separates  the  mass  of  the  earth  into  two  pans  -  one 
is  under  the  point  level  and  one  is  above  the  point  level.  Based  on  the  above  discussion  the 
disturbing  potential  generated  by  the  mass  under  the  point  level  can  be  represented  exactly  at  point 
P  by  an  expansion  similar  to  (4-9)  provided  the  downward  continued  potential  is  convergent.  The 
potential  of  the  mass  above  the  point  level  cannot  be  represented  exactly  by  an  expansion  like  (4- 
9). 


In  urder  to  simplify  this  problem,  we  take  the  planar  approximation.  Obviously,  the  potential 
of  the  mass  above  the  print  level  continued  downward  to  the  ellipsoid  is  given  by  (cf.  eq.  (4-2)); 
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where  Ah  is  defined  by  eq  (4-4);  V*  denotes  the  potential  on  the  ellipsoid. 
Upward  continue  the  to  the  point  level,  we  obtain 


(4-12) 

It  IS  the  same  as  eq  (  4  2)  At  the  same  time  the  potential  of  the  ma.ss  above  the  point  level  is  given 
by  eq.  (4  3i,  The  difference  between  the  Vp  and  Vp  i.s  given  by  eq.  (4-5).  Therefore,  the 
disturbing  potential  on  the  point  level  (which  is  partly  inside  the  earth)  represented  by  a  spherical 
harmonics  like  eq.  (4-9)  differs  the  true  disturbing  potential  by  a  quantity  oV  which  is  given  by  eq. 
(4-5)  Therefore  the  representation  of  the  disturbing  potential  by  an  expansion  of  spherical 
harmonics,  which  i.s  determined  by  using  the  gravity  anomalies  downward  continued  from  outside 
Earth's  surface  or  from  Earth's  surface  to  the  ellipsoid,  is  just  as  good  as  solving  the 
Molodensky's  problem  by  using  the  ADC  procedure. 

The  effect  of  the  topography  on  the  computation  of  the  geoid  from  a  spherical  harmonic 
expansion,  such  as  eq.  (4-9),  is  referenced  to  Wang  (1989). 
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5.  Conclusion 

The  e;mh'>;  snr*;-.ce  e;  „  coir.piic^tcu  suiface.  It  is  not  suitable  for  the  computations  of  the 
physical  geodesy.  Therefore  the  different  reference  surfaces  such  as  the  sphere,  the  ellipsoid  or  the 
point  level  have  been  used  for  solving  the  geodetic  boundai^'  value  problem.  It  is  clear,  if  such 
reference  surfaces  are  chosen,  the  analytical  downward  continuation  has  to  be  used  and  the  effect 
of  the  topography  has  to  be  considered.  Up  to  date  the  airborne  and  satellite  gradiometry  have 
been  reaching  great  improvement  and  will  supply  a  huge  gravity  gradient  data  for  the  geodetic 
boundary  value  problem.  To  process  such  data  the  reference  surface  has  to  be  chosen  and  the 
ADC  method  has  to  be  used.  Therefore  the  effect  of  the  topography  must  be  taken  into 
consideration. 

In  Chapter  2  the  topographic  correction  to  the  airborne  gravity  gradient  data  was  considered. 
In  this  chapter  the  Helmert’s  second  condensation  was  introduced  into  the  airborne  gradiometry. 
This  is  much  more  comple.x  than  it  is  in  the  classical  meaning  because  the  gravity  gradient  has  six 
components,  the  computations  of  the  topographic  correction  to  it  needs  much  more  computation 
effort.  It  was  found  that  the  Bouguer  plate  has  no  contribution  to  the  gravity  gradient.  The 
topographic  correction  is  mostly  due  to  the  terrain  of  the  earth. 

In  Chapter  3  we  introduced  the  remove- restore  technique.  Our  focus  was  still  on  how  the 
geoid  inside  the  earth  can  be  determined  by  using  the  method  of  the  analytical  downward 
continuation.  A  simple  way  to  do  this  is  to  obtain  the  disturbing  potential  at  the  flight  altitude  by 
processing  the  aerial  gradient  data,  and  then  downward  continue  it  to  the  sea  level.  Obviously, 
such  potential  is  not  the  original  one.  Some  correction  has  to  be  considered.  By  using  the  remove- 
restore  technique  we  found  that  the  topographic  correction  consists  mostly  by  a  term  -27tGpl^  and 
a  terrain  effect  which  can  be  very  small  after  the  smoothing.  The  first  term  is  important  for  the 
determination  of  the  geoid  in  a  mountain  area  and  the  terrain  effect  (eq.  (3-34))  is  very  rough  and 
becomes  small  when  some  kind  of  smoothing  are  applied.  Because  the  remove-restore  technique 
does  not  change  the  location  of  the  mass  of  the  topography,  the  indirect  effect  is  equal  to  zero.  In  a 
planar  approximation  the  change  of  the  potential  6V  and  its  Fourier  transformation  was  shown 
convergent  almost  everywhere,  except  at  the  infinite  point  of  the  circle  frequency  a)=«’.  This 
means  that  the  ADC  is  also  convergent  almost  everywhere.  This  is  imponant  for  developing  the 
geopotential  into  harmonic  series  up  to  high  or  very  high  degree  and  order. 

In  Chapter  4  the  topographic  effect  on  the  solution  of  Molodensky's  problem  by  using  the 
ADC  was  considered.  Theoretically,  the  topographic  effect  still  exists,  even  if  the  analytical 
downward  continuation  solution  is  equivalent  to  the  Molodensky's  solution  which  is  considered 
theoretically  perfect.  Because  the  topographic  effect  (eq.  (4-7))  can  be  very  small  after  smoothing, 
this  effect  can  be  neglected  in  the  practical  computations. 

The  same  topographic  effect  exists  also  for  the  geopotential  represented  by  a  spherical 
harmonic  expansion.  One  should  not  be  surprised  about  this  because  the  ADC  is  used  in  the 
solving  of  .Molodensky's  problem  and  the  determination  of  the  spherical  harmonics.  The 
discussion  in  Chapter  4  shows  that  a  simple  spherical  harmonic  function,  like  eq.  (4-9),  cannot 
represent  exactly  the  disturbing  potential  nearby  the  earth  and  on  the  eanh's  surface.  By  using  the 
ADC  method  we  can  obtain  a  spherical  harmonic  function,  it  represents  exactly  the  disturbing 
potential  on  and  above  the  boundary  sphere  of  the  earth,  provided  such  a  spherical  harmonics 
exists  (convergent);  it  represents  practically  well  the  disturbing  potential  nearby  the  earth  and  on  the 
eanh's  surface,  even  though  the  topographic  effect  (correction)  exists. 

In  practice  the  spnerical  harmonic  expansion  is  taken  to  a  finite  degree  and  order  and  some 
kind  of  smoothing  are  always  involved  in  the  practical  computation.  Therefore  the  topographic 
effect  can  be  very  small,  even  though  it  looks  very  rough  and  unstable.  Here  we  must  not  be 
confused  by  the  word  "topographic  effect".  The  topographic  effect  is  defined  as  the  difference  of 


rhe  true  geopotemial  on  the  point  level  and  the  potential  determined  bv  using  analytical  downward 
continuation.  Bectiuse  the  difference  is  caused  by  the  topography  of  the  earth,  we  cal!  this 
diffe’-'^ru'i*  topographic  effect  ( correction j.  'Ihis  effect  is  tot;ill\  different  from  the  classic 
topography-reduction  or  the  topographic  correction  to  the  spherical  harmonic  expansion  of  the 
geopotential  which  has  been  done  by  Sjdberg  (1988h 

Theoretically,  the  topographic  effect  is  imponant,  because  it  gives  the  idea  of  how  good  a 
spherical  harmonics  represents  the  disturbing  potential  on  and  outside  the  earth.  The  equivalent 
problem  is  how  good  can  the  disturbing  potential  on  the  eanh's  surface  be  determined  by  using  the 
ADC  method. 

There  is  only  one  gravity  field  but  it  can  only  be  approximated  by  using  the  numerically 
computed  different  gravity  models.  Based  on  the  results  in  Chapter  4  we  should  not  worry  about 
the  expanding  of  the  sphencal  harmonics  into  a  high  degree  and  order,  if  the  gravity  data  are  dense 
enough  and  accurate  enough  for  the  needs  In  such  a  case  the  sphencal  harmonic  expansion  of 
high  degree  and  order  will  represent  the  disturbing  potential  better  than  one  in  lower  degree  and 
order.  If  the  accuracy  of  the  spherical  harmonic  representation  is  so  high  that  the  topographic 
effect  has  to  be  taken  into  account,  the  topographic  effect  can  be  computed  with  a  proper 
smoothing. 

Note  that  the  analytically  downward  continuation  of  a  spherical  harmonics  into  the  eanh  is 
assumed.  Theoretically  it  is  difficult  to  prove  this  series  convergent  or  divergent  in  a  spherical  case 
even  though  it  was  shown  in  Chapter  3  that  the  ADC  is  convergent  almost  everywhere  in  a  plane 
approximation  case.  As  Moritz  (1980)  pointed  out,  such  a  problem  is  meaningless  in  practice.  We 
can  still  believe  that  a  spherical  harmonic  expansion,  like  eq.  (4-9),  can  represent  the  gravity  field 
as  well  as  we  wish.  The  derivation  in  Chapter  3  helped  us  to  heighten  our  confidence. 
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